Hydrogen transport in superionic system Rb3H(Se04)2: a revised cooperative 

migration mechanism 
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We performed density functional studies of electronic properties and mechanisms of hydrogen 
transport in Rb3H(Se04)2 crystal which represents technologically promising class M3H(X04)2 of 
proton conductors (M=Rb,Cs, NH4; X=S,Se). The electronic structure calculations show a decisive 
role of lattice dynamics in the process of proton migration. In the obtained revised mechanism of 
proton transport, the strong displacements of the vertex oxygens play a key role in the establishing 
the continuous hydrogen transport and in the achieving low activation energies of proton conduction 
which is in contrast to the standard two-stage Grotthuss mechanism of proton transport. Conse- 
quently, any realistic model description of proton transport should inevitably involve the interactions 
with the sublattice of the XO4 groups. 

PACS numbers: 71.15.Mb,71.20.-b,66.30.Dn 



I. INTRODUCTION 

Last years demonstrate a continuous increase of re- 
search activity in the field of hydrogen conductors. This 
can be explained by great technological perspectives of 
hydrogen conducting materials for applications in solid- 
state hydrogen fuel cells, hydrogen storage and elec- 
trochemical devicesi"— . A central problem in fuel-cell 
and hydrogen batteries technology is the development of 
cheap and efficient materials for electrochemical elements 
which can be used for the chemical storage and transfor- 
mation of chemical to electrical energy. In the functional 
properties of hydrogen devices, the transport of hydrogen 
through the conducting part to cathode plays a key role 
which makes the proton transport mechanisms one of the 
central problems of current experimental and theoretical 
studies^"—. 

In the wide range of proton conducting materials, the 
crystals with superionic phases M3lI(X04)2 (M=Rb,Cs, 
NH4; X=:S,Se) are of especial interest due to their high 
proton conductivity coefficients which in the high tem- 
perature disordered state can reach values of the order 
10~^ — 10""'^ Q.~^cm~^. In these phases, the quasi-two- 
dimensional proton transport is characterized by a dy- 
namically disordered network of virtual hydrogen bonds 
dynamically establishing between the oxygens 0(2) of 
nearest groups XO4. Although on the first sight the 
transport of protons in hexagonal (001) planes seems to 
be quite a strainforward process explained in terms of the 
two-stage Grotthuss mechanism^^^'^, the real mechanism 
of proton migration is still not suflBciently well exam- 
ined and understood. The reason for this difficulty is a 
complex cooperative character of proton migration which 
should strongly interact with the dynamics of the ionic 
environment, the fact supported by numerous experi- 
mental Raman, Infrared, x-ray and NMR spectroscopy 
studiesii"— . 

In the theoretical worksi^ii^"— , the two-stage Grot- 



thuss mechanism is frequently considered within phe- 
nomenological free-particle approches which typically 
contain two stages: (1) translation of proton within a 
hydrogen bond between the nearest ionic groups XO4; 
(2) inter-bond proton transfer related to reorientations of 
the ionic groups XO4. Howerer, the experimental stud- 
ies demonstrate that the transport of protons is a highly 
complex process which can be characterized as a transfer 
of hydrogen ion in the course of the reactions of cre- 
ation and breaking of IIXO4 bonds. As a consequence of 
such reacting environment with dynamically redistribut- 
ing electronic density, the protons in conducting materi- 
als cannot be described as almost free particles. 

Despite the existing theoretical studies of the low- 
temperature proton ordered phases and proton conduc- 
tion in M3H(X04)2 by phenomenological modelai^"— , 
up to date there is no first-principle-investigations of 
these systems. The recent first-principle molecular dy- 
namics studies of superprotonic phase transition in a 
related superionic system CsHS04^2, were focused on 
the structural transformations and reorientational mo- 
tion ofHSO^ groups, without the studies of the proton 
migration mechanism. Another ab-initio study of the 
proton conduction in CSHSO4 has beed devoted to the 
estimation of diffusion barriers, without a detailed con- 
sideration of the role of oxygen network in the proton 
transport mechanism^!. The lack of the first-principle in- 
vestigations of the transport mechanism in M3H(X04)2 
can be explained by highly complicated crystal struc- 
ture of these systems involving several structural trans- 
formations which makes ab-initio calculations extremely 
time-consuming and challenging. Our present studies of 
the hydrogen transport are based on the density func- 
tional theory (DFT) with the use of the linearized aug- 
mented plane wave method (LAPW) implemented in the 
WIEN2k packagers. The LDA (Local Density with non- 
linear core corrections) and GGA (Generalized Gradient) 
approximations of the DFT developed by Perdew, Burke 
and Ernzerhof^^, and Perdew and Zunger^^ have been 
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employed in our calculations. In the GGA, the additional 
gradient terms of the electron density are added to the 
exchange-correlation energy and its corresponding po- 
tential, which is beyond the local density approximation 
(LDA) of the DFT. In our studies of the proton migra- 
tion paths, we also employ a nudged-elastic-band method 
within the pseudopotential approach implemented in the 
Quantum Espresso package^^ which allows to perform the 
dynamical relaxation of the atomic surrounding along the 
steps of the proton conduction processes. 

As an attempt to account for the coupling of mobile 
protons with the dynamical distortions of the groups 
XO4, the existing microscopic description of the proton 
transport considers the interaction of protons with the 
displacements of the vertex oxygens 0(2) of XO4 involved 
in the hydrogen bonding^-. Such a coupling is presented 
in terms of the interaction of proton charge with the 
optical phonon modes corresponding to the anti-phase 
stretching vibrations of the 0(2)-ions and is classified 
in terms of "protonic polaron" concept. The calculated 
in this approach coefficients of the proton conductivity 
are consistent with the experimental measurements per- 
formed in the superionic phases^^. However, the pola- 
ronic concept introduced in Ref. [2l| involves a number 
of phenomenological parameters like polaron binding en- 
ergy which have been estimated by the fitting of the cal- 
culated conductivity to the experimetally measured val- 
ues. To verify the role of the dynamical distortions in the 
proton migration mechanism, in this work we present a 
detailed DFT study of the mechanism of hydrogen mi- 
gration in the system Rb3H(Se04)2, which belongs to 
the M3H(X04)2-class. 

To investigate the role of the proton environment, two 
aspects of the hydrogen transport have been analyzed. 
The first aspect is related to the structure and position- 
ing of hydrogen-bonded oxygen ions 0(2) in the dynami- 
cal network. The present DFT calculations demonstrate 
the strong displacements of 0(2) from the high-symmetry 
positions on the three-fild axis in the equilibrium state, 
which allows to verify the concepts of the 0(2)-structural 
positioning discussed in the previous x-ray and NMR 
studiesi^i^S. As the second aspect, we have analyzed 
the migration paths of protons in the configuration of 
instantaneously relaxed atoms. The performed DFT- 
calculations allow to derive the energy profiles along the 
proton migration paths and to investigate the displace- 
ments of each atom involved in the proton surrounding. 
As a result, the present studies give a possibility to obtain 
the first-principle estimates of the hydrogen-bond lengths 
and energetic parameters involved in the phenomenologi- 
cal modeling and to verify the concept of protonic polaron 
introduced in the previous phenomenological description 
of proton transport. Our findings also allow to improve 
several experimental conclusions obtained in the previous 
x-ray and NMR experiments with Rb3lI(Se04)2^>^, in 
particular to clalify the positioning of the vertex oxygens 
0(2) during the proton migration and to examine the 
proton positions on the hydrogen bond at the intrabond 



proton transfer. As a result, the present DFT results 
demonstrate a central role of oxygen displacements for 
the migration of hydrogen and emphasize in this way a 
key role of the proton environment for the establishing 
the proton transport. 

The rest of the paper is organized in three sections. In 
Section|TTl the theoretical methods employed in our stud- 
ies are discussed. Section IIIII presents the results of the 
electronic structure calculations and of the calculations 
of the proton migration paths and oxygen positioning. In 
Section IIVI we summarize the main results and conclu- 
sions of our work. 



II. THEORETICAL MODELS AND METHODS 

A. Structural concepts of hydrogen-bonded 
network 

In the high-temperature disordered phase of 
Rb3H(Se04)2 schematically shown in Fig. [IJa), hy- 
drogen transport occurs in a network of "virtual" 
(equally probable) hydrogen bonds which are occupied 
by a proton with equal probability 1/3. This network 
(see Fig. [Ijb) for details) is established between the 
vertex oxygens 0(2) of the tetrahedra Se04 which can 
be chemically bonded to the three nearest tetrahedra by 
hydrogen bonds. 

The experimentally deduced structure of the network 
of hydrogen-bonded tetrahedra at T = 456 K projected 
on (001) plane is represented by the electronic difference 
density map in Fig. [2] This map has been calculated 
from the x-ray diffraction data using SHELXL-93 pro- 
grams. The central feature detected in Fig. [2] is a disor- 
der of the hydrogen-bonded vertex oxygens 0(2) (shown 
by red dots) between three structurally equivalent posi- 
tions. The hydrogen bonds are schematically presented 
by the solid lines where the double black dots denote the 
positions for proton in a double- well potential. 

In the earlier work by Baranov et alj^^, a different con- 
cept of the structure of hydrogen-bonded network is pro- 
posed for the analysis of the x-ray diffraction data. In 
this concept, the oxygens 0(2) are statically located in 
high-symmetry positions on three-fold axis, which also 
implies a static character of the of Se04 tetrahedra and 
leads to a possibility of the movement of free protons in 
the disordered superionic phase. Therefore, the arrange- 
ment of the vertex oxygens 0(2) is still an open question 
which has a prime importance for the understanding the 
nature of proton transfer mechanism in these systems. 

In the first concept of the structural arrangement of 
0(2) in the superionic phase of Rb3H(Se04)2 signified 
in Ref. JLSj, the oxygens 0(2) are located in the high- 
symmetry positions on the three-fold axis and form weak 
hydrogen bonds with the length of the covalent bond 
0(2)-II approaching 1.76 A. In distinction from the high- 
symmetry case, the second type of oxygen arrangement 
discussed in Ref. [TsIITtI assumes the distortions of 0(2) 
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FIG. 1: (a) Schematical presentation of the crystal structure 
of Rb3H(Se04)2 in the superionic phase in projection on the 
hexagonal (001) plane, where the rhombohedral unit cell of 
R3m symmetry is doubled along the (lOO)-direction. One of 
the virtual hydrogen bonds with the proton centered on the 
bond is shown for simplicity. Here the small orange circles 
denote the protons, the oxygen atoms are indicated by blue, 
Rb atoms by green and Se by brown circles, (b) Atomic con- 
figuration of hydrogen-bonded network of Se04 tetrahedra in 
superionic phase in the proton conducting layer correspond- 
ing to a slab (—0.13c < z < 0.13c) in the hexagonal (x,y) 
plane, with the lattice constant c = 22.629A. The orange cir- 
cles indicate the possible positions for the disordered protons 
occupied with a probabibity 1/3. 



from the three-fold axis towards the hydrogen at the for- 
mation of the 0(2)-H-0(2)' bond. To prevent a break of 
high symmetry in the disordered superionic phase caused 
by such distortions, a dynamical disorder of each 0(2) 
between three equally probable positions is assumed. 
These positions are related to each other by a rotation by 
120 deg around the three- fold axis (see Fig. [2] where the 
three equivalent positions of distorted 0(2) are marked 
by red dots)). 

To analyze the equilibrium positions of the oxygens 
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FIG. 2: Electron charge density maps of the nearest neigh- 
bouring Se04 groups of Rb3lI(Se04)2 in the disordered supe- 
rionic state at T = 456 K. Here the red dots show the possible 
positions of the vertex oxygens 0(2) around the three- fold axis 
and the black dots correspond to the location of disordered 
protons on 0(2)-H-0(2)' hydrogen bonds. 



0(2), we performed the GGA calculations of the total 
energy of the system for the different values of the dis- 
placements A of the oxygen 0(2) from the positions of 
the three-fold axis towards a migrating proton, a process 
which is schematically shown in Fig. [TJb). In the cal- 
culations of the electronic structure of Rb3H(Se04)2, the 
lattice constants of the unit cell of the trigonal R3m sym- 
metry in hexagonal coordinates were fixed to the values 
obtained from the structural analysisi^: a — 6.118 A 
and c — 22.629 A. Technical details include a GGA 
method developed by Perdew, Burke and Ernzerho^ 
on a 4 X 2 X 3 fc-point grid representing a mesh of 
40 fc-points in the first Brillouin zone. The compari- 
son of the total energy calculated for the mesh of 80 fc- 
points gives the difference between the total energies 
Etot{80) - Etoti^O) = 0.00136 eV which allows us to use 
the 40 fc-point mesh for the GGA-calculations. 

The full-potential method^^ employs an expansion of 
the electronic potential inside the atomic spheres via a 
full number of the electronically occupied local orbitals 
(Rb, Se, O and H atoms in our case), whereas the po- 
tential outside the spheres is constructed as a plane- wave 
expansion with a plane- wave cutoff given by 667 eV. Such 
a method allows to achieve high accuracy in the descrip- 
tion of the electronic structure, total energy and elec- 
tronic density of states. The dynamical disorder in the 
superionic state cannot be directly described in the static 
density functional calculations. As a consequence, in a 
generated model unit cell, only a single hydrogen bond is 
selected from the three virtual hydrogen bonds near each 
Se04 group which leads to a lowering of the symmetry of 
the system from the rhombohedral to monoclinic. In this 
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case, the generated structure can be considered as a static 
snapshot of the dynamically disordered superionic state. 
In the studies of the hydrogen transport mechanism we 
also considered a doubled unit cell, where the initial ori- 
entations of the two hydrogen bonds in the cell are chosen 
randomly, each bond can be dynamically broked and cre- 
ated in a chain of intermediate configurations generated 
by the nudged elastic band (NEB) niethod^®''^^. 

In the studies of the transport mechanism, each atomic 
configuration has been structurally relaxed. The relax- 
ation involves the optimization of the atomic positions 
by the minimization of the total energy and forces acting 
on the atoms in the generated unit cell. For the struc- 
tural optimization, we use a Newton scheme described in 
Rcf. 31. A comparison of different relaxed static config- 
urations allows to calculate the energy barriers for the 
proton migration. With the obtained zero-temperature 
results, we also analyze the influence of temperature on 
the activation energies of the proton transport. 



B. Modeling the hydrogen migration paths 



properties calculated by these two approaches are consis- 
tent, we conclude that the PZ-pseudopotential approach 
implemented in the QE-code is sufficient to capture the 
physical features of the hydrogen migration in our sys- 
tem. 
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To analyse in details the paths for the hydrogen 
migration in the relaxed surrounding, we employed 
the nudged-elastic-band (NEB) method implemented in 
the Quantum-Espresso (QE) Package of the DFT cal- 
culations with the use of plane-wave basis sets and 
pseudopotentials^Si^. In these calculations, for the 
atomic cores of Rb, H, Se and O we employ the 
Perdew-Zunger (PZ) norm-conserving pseudopotentials 
with nonlinear core corrections^!. Our choice of pseu- 
dopotentials is justified by the facts that all characteris- 
tic features of the electronic structure calculated with the 
use of the full potential Wien2k GGA code (Fig. ^a.)) 
and obtained with the PZ-pseudopotentials (Fig. |3jb)) 
are consistent in the energy window (Ep — 7.0eV; Ep + 
leV), where Ep is the Fermi energy. Specifically, both 
approaches allow to describe the characteristic DOS- 
peaks related to the 0(2)-2p orbitals (the energy range 
between Ep and Ep — 5 eV in Fig. [3]) and the DOS-peak 
related to the Is-orbital of H at the energy value about 
—6.1 eV (Fig. |3l[a)) which corresponds to the peak at 
about —6.2 eV in Fig. [3Kb). In the GGA-approach, the 
energy gap between the highest filled and lowest empty 
electronic states is about 3.5 eV (Fig.[3Ka)) which is close 
to experimentally discussed gap about 4 eV— in related 
compounds. The GGA gap is larger then the gap about 
1.3 eV obtained by the LDA-method (Fig. in^b)) due to 
the usual underestimation of the gap in the LDA cal- 
culations, which does not influence the results of the 
NEB studies based on the calculations of the total energy. 
Also, we would like to note that a detailed comparison 
of the electronic and optical properties of several com- 
pounds like Rb20 and Li20 related to the constituents 
of Rb3lI(Se04)2 calculated in the LDA and GGA ap- 
proaches shows a good consistency with the available ex- 
perimental results'^'^. As the low-temperature electronic 



FIG. 3: Total density of states calculated (a) in the GGA 
(Perdew, Burke, Ernzerhof, Ref. (2^ ) using the full potential 
approach of Ref. [2^ and (b) in the LDA using norm-conserving 
atomic pseudopotential approach with nonlinear core correc- 
tions implemented by Perdew, Zunger, see Ref. [27l . The zero 
energies correspond to the Fermi levels. 

In contrast to the almost free motion of the hydrogen 
implemented in a wide range of the models of hydrogen 
transport, in the present work we consider the migration 
of hydrogen as a cooperative process, the feature which 
can be described by the NEB-approach22i2i. Specifically, 
for each stage of the hydrogen transport, this process 
involves a relaxation of the atomic positions and of the 
distances between the different atoms. In our analysis, to 
better relax the positions of the atoms neighbouring to H, 
we perform a doubling of the unit cell which now includes 
four Se04 groups and two protons (see Fig. [DJa) for de- 
tails) . For each position of the hydrogen in the unit cell, 
the coordinates of all neighbouring atoms were relaxed 
until the forces acting on the atoms reached their min- 
ima. In these calculations, we use the plane-wave cutoff 
1020 eV and the energy cutoff for charge and potential 
given by 2040 eV. 

In the NEB-approach, the relaxation of the atomic po- 
sitions along the hydrogen migration path is performed 
by the minimization of the total energy of each interme- 
diate configuration (image). These images correspond to 
different positions of H on the migration path and pro- 
duced by the optimization of a specially generated object 
functional (action) with the consequent minimization of 
the spring forces perpendicular to the path. In our calcu- 
lations, the convergence criteria for the norm of the force 
orthogonal to the path is achieved at the values below 
0.05eV/A. 
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III. RESULTS AND DISCUSSION 



A. Vertex oxygens 0(2) 



The experiments (see Ref. |29I ) show that in 
M3H(X04)3-crastals the temperature of the superionic 
phase transition Ts does not depend on the deuteriza- 
tion. Moreover, the neutron diffraction resuhs^ do not 
distinguish two distinct hydrogen positions on the bonds 
and justify a representation oh hydrogen atom in flat 
single-minimum potentials. Consequently, in the present 
studies of the positioning of 0(2), the initial proton po- 
sitions are fixed to the centers of hydrogen bonds. In 
Fig. [1] such a single- well approximation is identified by 
the centering of proton between the oxygens 0(2). 

To study the modification of the total energy caused 
by the relaxation of the vertex oxygens, the positions of 
0(2) on a selected hydrogen bond were shifted by a dis- 
tance A, while the coordinates of other atoms were kept 
in the undistorted high-symmetry state. Fig. |4] demon- 
strates the variation of the total energy of the system 
E{A) - E{A = 0) with the increase of the 0(2) dis- 
placement parameter A. One can see that the lowest 
value of E{A) corresponds to A = 0.5 A which im- 
plies a stabilization of a short hydrogen bond 0(2)-H- 
0(2)' of a length 2.54 A where the length of the cova- 
lent bond 0(2)-H equals 1.27 A. It is remarkable that 
despite the central-point approximation for the proton 
within the hydrogen bonds, the obtained length of the 
H-bond is typical for double-well H-bonds and is con- 
sistent with the results of x-ray and neutron scattering 
studiesi^iiii^, which becomes clear from the comparison 
of different data with the present results presented in 
Table ID Consequently, our findings clearly support the 
second concept of the oxygen-distorted network structure 
which is characterized by the two types of dynamical dis- 
order in the superionic phase: (i) the statistical disorder 
of the hydrogens between the virtual hydrogen bonds in 
hexagonal planes and (ii) the disorder of the vertex oxy- 
gens between three symmetry-related positions around 
the three- fold axis indicated in Fig. [IJb). 

To demonstrate the influence of the displacements of 
0(2) on the electronic properties, we have analyzed the 
projected atomic densities of electronic states for A = 
and A = 0.5 A which are presented in Fig.[S] The consid- 
ered two types of the structural positioning of 0(2) can 
be characterized by a substantial electronic hybridization 
of H Is orbitals and 2p orbitals of 0(2). In Fig. [SJ the 
s-p hybridization can be identified by the density peaks 
at the energies about —3.5 eV and in the range (—5.5)- 
(—6) eV below the Fermi energy (for A — 0.5 A, these 
peaks are marked by a dashed line). The DOS peaks at 
about —18 eV imply also a strong hybridization between 
H Is and 0(2) 2s orbitals. The displacements of 0(2) 
from the high-symmetry positions lead to a shortening 
of the 0(2)-H distances and consequently to a signifi- 
cant increase of the s-p hybridization peaks at —6 eV. 
Moreover, the 0(2)-distortions result in the emergence 
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FIG. 4: Total system energy as a function of the 0(2) distor- 
tion parameter A. Here A is the displacement of 0(2) with 
respect to its high-symmetry position on the three-fold axis. 



of new strongly hybridized state at —9.5 eV. In addition, 
the displacements of 0(2) induce the distortions of the 
Se04 tetrahedra which results in the lowering of the lo- 
cal symmetry and in the split of high hybridized peaks 
around —18 eV (Fig. [51 case A — 0.5 A). In the con- 
sidered model unit cell, instead of the virtual H-bonds 
we consider a number of selected configurations with the 
fixed H-bonds near each Se04 group. In a real statisti- 
cally disordered system one may expect that the obtained 
split of the hybridization peaks will be smoothed by the 
dynamical disorder. 




Energy (eV) 

FIG. 5: Projected density of states of the oxygens 0(2) and 
hydrogens calculated by the GGA method. The cases A = 
and A — 0.5 A correspond to the hydrogen bonds between 
Se04 groups (a) with high-symmetry positions of 0(2) on the 
three-fold axis and (b) with the displacements of 0(2) from 
the three-fold positions towards the proton. The zero energy 
corresponds to the Fermi level. 

The redistribution of the electronic charge density due 
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to the distortions of 0(2) causes the accumulation of the 
electron charge on the obtained short hydrogen bonds, 
the effect observed in the contours of the difference charge 
density in Fig. [51 Consequently, in the superionionic 
phase the formation of the strong virtual hydrogen bonds 
should be also identified by a remarkable increase of the 
electron density between the 0(2) ions, which is consis- 
tent with the experimentally obtained electron density 
map of Fig [T] 



Aq=0 a ^o~°-^ ^ 




coord i nate X (A) coord i nate x (A) 

FIG. 6: Electronic difference density contours (crystalline 
renormalized by the superimposed atomic densities) for the 
hydrogen bonds 0(2)-H-0(2). Here the oxygens 0(2) are lo- 
cated on the three- fold axis (left panel), and displaced from 
this axis towards the hydrogen on the bond (right panel). 



B. Rotational stage in the full potential 
calculations 

Based on the results obtained for the positions of the 
hydrogen-bonded 0(2), it is possible to study the mecha- 
nisms and paths of the hydrogen migration. In our anal- 
ysis, the special focus is put on the rotational stage of 
the Grotthuss transport mechanism. In its classical de- 
scription, this stage involves the rotational motion of the 
HSe04 group around the three-fold axis determined by a 
rotation angle a (see Fig. [7l[b) where the rotation stage 



is schematically shown), and the consequent formation of 
a new hydrogen bond between the rotated group and the 
nearest Se04 tetrahedron approached by this rotation. 

First, we analyze a possibility of the rotational step 
as an independent stage of the proton transport mech- 
anism. In this case, the hydrogen with the covalently 
bonded oxygen 0(2) is allowed to rotate around the near- 
est SeOa-group while the positions of all the surrounding 
atomic groups are fixed. In such a rotational stage, we 
have calculated the variation of the total energy of the 
system due to the rotations. Similar to the studies of the 
0(2) network in the previous section, we consider here 
two possible structural concepts for the rotational mo- 
tion. In the first high-symmetry case, the oxygen 0(2)^ 
of the II(Se04)"'^ group is fixed on the three-fold axis 
(A = 0), whereas in the second configuration the oxy- 
gen 0(2)^ is relaxed (displaced) by a distance A = 0.5 A 
from the three-fold axis and rotates jointly with the hy- 
drogen by the angle a (the rotation of 0(2)^-11 covalent 
bond) . The resulting dependences of the total energy on 
a are presented in Fig. [T] 

It is remarkable that both, high-symmetry and 0(2)^- 
H-rotated configurations are characterized by high- 
energy barriers (about 2.6 — 3 eV) between the initial po- 
sition of proton and the transition state (TS) (Fig.lTJa)) 
which corresponds to the intermediate value of the rota- 
tion angle a about 7r/3. The high energy barriers ahso- 
lutely exclude any independent motion of proton in the 
conducting system, a concept which is widely applied for 
the free electrons in metals. We also note the obtained 
values for the energy barriers are af the same order as the 
rotational barrier 1.54 eV calculated in Ref. f23 for the 
system CsHSe04 which prompts on the common physi- 
cal aspects of the transport mechanisms in these classes 
of proton conductors. 

Despite the high energy of the TS state, the dashed 
energy profile in Fig. [7] demonstrate a possibility for a 
partial rotation of the covalent bond 0(2)-'^-H by a = 
20 — 25 deg which has considerably low barriers of about 
0.3 — 0.7 eV. In the superionic state, such a partial rota- 
tion can be also accompanied by the displacements of the 
oxygen 0(2)^ of a neighbouring group (Se04)^ towards 
the hydrogen which would lower the energy barrier and 
could result in the formation of the new hydrogen bond 
0(2)^-H-0(2)^ with the consequent proton transfer to 
the group (Se04)^. As a conclusion, the rotational mo- 
tion of protons should also involve the relaxation of the 
oxygens 0(2) of more distant Se04 groups. 



TABLE I: Hydrogen bond length parameters (in A) in 
the superionic phase of Rb3H(Se04)2 determined from the 
experiments and calculated by the 0(2)-relaxation in the 
full-potential GGA. 



GGA (0(2)-rolaxation) Ref.J^ Ref. 30 Ref.il 
2.54 2.51 2.594 2.58 



C. Hydrogen migration as a cooperative process 

In the studies of hydrogen migration in a fully relaxed 
surrounding, we employ the NEB-mathod discussed in 
section IIIBI Fig. |8] shows the snapshots of the migra- 
tion of the hydrogen HI between the complexes (Se04)^- 
(Se04)^ and (Se04)^-(Se04)^ calculated by the NEB 
method. The obtained proton migration path can be 
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FIG. 7: (a) Energy profiles for tfie rotation of liydrogen 
around the group (Se04)^. The case A = corresponds to 
the oxygen 0(2)^ located on the three-fold axis and fixed dur- 
ing the rotation of H. For A = 0.5 A, 0(2)^ is displaced by 
a distance A and rotates jointly with the hydrogen, (b) She- 
matical view of the rotation stage of the transport mechanism. 
The short arrows indicate the displacements of the oxygens 
0(2) from the three-fold positions towards the proton, the 
large circles denote SeOa groups, a is the rotation angle of 
0(2)-H covalent bond. 



described via several steps: 

(1) the 0(2)^-111 covalent bond is rotated around 
SeOj-group which corresponds to the variation of the ro- 
tation angle a from its initial value 7r/3 to the value a — 
(part (a) in Fig. |8]is schematically shown in Fig. IH^c)). 
We note that this rotation fully corresponds to the rota- 
tional motion of H studied in the previous section IIII Bl 
within the full potential approach, whereas the rotation 
in section IIII Bl is performed in the opposite direction 
shown in Fig. with the angle a changing from 

to 7r/3. The comparison of the corresponding energy 
profiles calculated by the NEB-method ('Fig. Elja), HI 
coordinate r in the range < r < 0.8A) and in the 
full-potential approach (Fig. [S]) shows that despite the 
similar parabolic form, the energy barrier in the QE-case 
AE{QE) = £'(0) - E{r = 0.8) is about 1.0 eV which is 
much lower then the full-potential GGA-value obtained 
from the difference AE{FP) = E{a = tt/S) - E{a = 
0) = 2.72 eV. This comparison clearly demonstrates the 
key role of the structural relaxation in the NEB-approach 



for the lowering the energy barriers for the proton trans- 
port. 




V« Hlf (1) * 

. • (a) 



(2) 
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« r^* (c) 



T ^ 
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FIG. 8: Snapshots demonstrating hydrogen migration be- 
tween the neighbouring groups (Se04)^ and (Se04)'^. Here 
Se, O and H atoms are represented by yellow, red and blue 
circles. Sections (a) and (b) show configurations with rotated 
(Se04)^; (c) resulted from the transfer from (1) to (2) and 
rotation around (Se04)^, configuration (d) is obtained after 
the transfer of hydrogen between two tetrahedra (Se04)^ and 
(Se04)^ 



In the further rotation of 0(2)^-111 the HI coordi- 
nate r changes from r — 0.8 to r = 2A which leads 
to a small increase of the total energy by about 0.2eV in 
Fig. Ufa). In Fig. EJa), the minimum of the total energy 
at r = 0.8A corresponds to the nearly zigzag alignment of 
the 0(2)-H-dimers (configuration presented in part (b) of 
Fig. similarly to the low-temperature ordered phase 
HI stabilized below the superionic transition tempera- 
ture, see Refs. [l5ll20ll33 for the detailed description of the 
structure of this phase. The rotation part in the range 
0.8 < r < 2 A corresponds to the rotation stage of the 
Grotthuss transport mechanism. 
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FIG. 9: Energy profiles for a migration of iiydrogen between 
tlie tetraiiedra (Se04)^ and (Se04)^. Here tlie part (a) corre- 
sponds to tlie rotation of the covalent bond Hl-0(2) around 
the group (Se04)^ and the part (b) corresponds to the transfer 
of HI to the group (Se04)'^. (c) Schematic view of the migra- 
tion path of HI in the projection on the hexagonal (OOl)-plane 
where the small symbols (a)-(d) correspond to the images in 

Fig.m 



In Fig.[9l^a), the energy barrier for the hydrogen rota- 
tion step approaches 0.2 eV. Here the rotation barrier is 
estimated from the energy dependence between r — 0.8 
(proton position indicated by the dashed line in Fig.[3Jc)) 
and r = 2A. 

(2) The hydrogen is transferred from the oxygen ion 
0(2)1 (pg^j.^ ^1^) Fig.© to the oxygen 0(2)^ (part (c) 
in Fig.[5|). Fig. [TU] presents the displacements of the oxy- 
gens 0(2) and the changes of the 0(2)-0(2) distances 
during the migration of the hydrogen derived from the 
DFT(NEB)-calculations. The proton transfer stage oc- 
curs in the range (2 < r < 3) of the proton coordinate. It 
is remarkable that this range of r is characterized by sig- 
nificant displacements of 0(2)^ and 0(2)^ by 0.2-0.25 A 
towards the hydrogen and by the corresponding strong 
decrease of the distance between 0(2)^ and 0(2)^ from 
its inital value about SA to the value 2.5lA, which is very 
close to the results of the full-potential GGA studies (see 
Table IJ). In Fig.[TUl the occurrence of the transfer of the 
hydrogen between 0(2)^ and 0(2)^ is identified by the 
minimal 0(2)-0(2) distance and maximal displacement 




0.4 0.8 1.2 1.6 2 2.4 2.8 3.2 3.6 4 
HI path coordinate r (A) 




HI path coordinate r (A) 

FIG. 10: (a) Displacements of the vertex oxygens 0(2) and (b) 
the change of the distance between the nearest 0(2) along the 
migration path of HI. The red arrows identify the transition 
state of the proton between the hydrogen atoms 0(2)^ and 
0(2)2. 



of 0(2)1 which is indicated by the red arrows in the plots. 
Table HIl presents the maximal (rfocf ) minimal (rfoo ) 
lengths of the hydrogen bond during the transfer step cal- 
culated in the NEB approach which are compared to the 
corresponding experimentally obtained values. Despite 
the good agreement with the experimental maximal bond 
length, the calculated value of dgo = 2.51 A is in coin- 
cidence with the estimates of Ref. [l^ for the superionic 
Rb3H(Se04)2 but is larger then the value 2.41 A dis- 
cussed in Ref. [13 for the similar system (NH4)3H(Se04)2. 
The difference is possibly related to the influence of the 
dynamics of NH^" groups in the superionic state which 
leads to additional softening of the structure and to the 
consequent increase of the 0(2) displacements 

(3) The 0(2)^-Hl covalent bond rotates and the hydro- 
gen HI moves towards the group (Se04)^. The process 
of the 0(2)^-Hl rotation in the group (Se04)^ is accom- 
panied by the elongation of the covalent bond between 
the oxygen 0(2)'^ and another proton H2 (shown in the 
left corner of part (c)) and by a strong displacement of 



TABLE II: Hydrogen bond parameters (minimal and maximal 
0(2)-0(2) distances dgg' and docT in A), protonic polaron 
transport parameter W (eV/A) and the transport activa- 
tion energy iJa(eV) in the superionic phase of Rb3H(Se04)2 
calculated by the NEB-method and estimated experimentally. 



method 


imin 
"OO 


^max 
"DO 






NEB (full relaxation) 


2.51 


2.7 


2.55 


0.33 


Refs. 


2.4 


2.7 




0.37 


Ref. 


2.51 


2.67 




0.49 
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the oxygen 0(2)^ towards the hydrogen HI (part (d) in 
Fig.H]). The displacement of the oxygen 0(2)'^ leads to a 
decrease of the distance between 0(2)^ and 0(2)'^ from 
its initial value 4.2 A to the value 3.22 A (Fig.fTUJb)) and 
is essential for the successful hydrogen migration to the 
group (8604)"^. At the transfer of HI between 0(2)^ and 
0(2)'^, the lengths of the bonds Se-0(2)'^ (vertex oxygen) 
and Se-O(l)'^ (basis oxygen) in the group (SeO)'^ are also 
changed by about 0.08-0. 16A. After the proton is trans- 
ferred from 0(2)2 q(2)3^ ^^e length of the Se-0(2)2 
bond is returned to its previus value 1.76A. 




FIG. 11: Electronic density contours showing the redistri- 
bution of the valence charge density of the hydrogen bond 
0(2)^-H-0(2)^ during the intrabond proton transfer (a)-(d) 
and during the the rotational motion of H covalently bonded 
to 0(2)^ (d)-(f). The valence charge has been obtained by the 
integration of DOS in the energy window {Ef — 4 eV; Ef) 
where Ep is the Fermi level. 

Fig. [TT] demonstrates the redistribution of the valence 
charge density during the migration of the hydrogen. The 
intrabond proton transfer step between the vertex oxy- 
gens 0(2)^ and 0(2)^ corresponds to the contours (a)-(d) 
of the plot. One can see that the hydrogen transfer has a 
character of a continuous redistribution of the electronic 
charge density between the two different positions near 
0(2)^ and near 0(2)^ detected in the plots (b) and (c). 
This finding allows to improve the conclusions of the neu- 
tron scattering studies in Ref. 'so' about the absence of the 
two distinct proton positions within a hydrogen bridge 
in the conducting phase of Rb3H(Se04)2 which can be 
explained by the limited resolution of the experiments. 
Furthermore, in the rotational step (plots (d)-(f)), the 
rotational displacements of the oxygen 0(2)^ follow the 
proton rotation and lead to significant modifications of 
the charge density contours near 0(2)^ which obtain a 
distinct dumbbell character. 

As the consequence of our findings, we propose to con- 



sider the obtained in our studies significant displacements 
of the oxygens 0(2) during the migration of the hydrogen 
as an additional feature of the proton transport mecha- 
nism. It is also remarkable that the well established in 
the literature standard Grotthuss proton transport mech- 
anism contains two main (rotational and translational) 
stages and does not involve the oxygens displacements as 
a decisive factor in the proton migration process, which 
is in contrast to the results of the present work. 

To see the role of 0(2) displacements in the establish- 
ing the hydrogen transport, we examine more closely the 
change of the total system energy along the proton mi- 
gration path presented in Fig. 13 

Fig. [S] shows that the energy barrier for the hydro- 
gen rotation step approaches 0.2 eV and is substantially 
lower than the barrier for the proton transfer between 
neighbouring groups which is about 0.6 eV in Fig. [HJb). 
The relation between the obtained rotation and transla- 
tion barriers is consistent with the conclusion about the 
decisive role of the intrabond proton transfer in the de- 
termining the activation energy for proton conductivity 
reported in Refs. lolfTTl 

In M3H(X04)2 systems, the dominant role of the in- 
trabond translational step in the formation of the trans- 
port activation energy is in contrast to the quasi-one- 
dimcnsional proton conductors MHXO4, where the pro- 
ton conduction process is determined by the rotational 
energy barriers. This fact is also supported by the re- 
cent DFT studies of CSHSO42"' where the intrabond pro- 
ton barrier 0.16 eV was found to be substantially lower 
than the rotational barrier 0.52 eV, and is also consis- 
tent with the ab-initio molecular dynamics calculations 
of the related systems which show strong fluctuations 
of the hydrogen bond lengths between 2.4 and 2.8 A 
during the proton migration and fast almost barrierless 
proton transfer on the hydrogen bonds characterized by 
the extremely short lengths about 2.4 A. The fast intra- 
bond proton transfer in the quasi-one-dimentional pro- 
ton conductors is similar to the proton transport proper- 
ties in water mixtures where the recent ab-initio molecu- 
lar dynamics simulations demonstrate the existence of 
two different time scales and identify the faster pro- 
cess with the intrabond proton transfer, while the slower 
time scale corresponds to local rearrangements of the 
hydrogen-bond network—. In this context, the unique 
two-dimensional character of the transport mechanism 
in the M3H(X04)2 introduces qualitative changes in the 
translational and rotational energy barriers with a con- 
sequent importance of the intrabond transfer in the de- 
termination of the transport energy barriers. 

The change of the total energy profile along the trans- 
fer path in Fig. IHJb) can be directly compared with the 
corresponding displacements of hydrogen bonded 0(2) 
in Fig. [TUT a) . From this comparison, we obtain that the 
increase of the total energy up to — 0.55 eV is accom- 
panied by the 0(2)-displacements Aq = 0.25 A. In the 
phenomenological description of the proton transport^i, 
the anti- phase displacements of 0(2) induced by the 
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hydrogen-bonded proton have been classified in terms 
of the "proton-polaronic" effect where the protonic po- 
laron is formed due to the interactions of proton with the 
anti-phase stretching displacements of the nearest 0(2)- 
oxygens and the energy for the formation of the protonic 
polaroni?o = (WF)^/2A/(cjo)^ depends on the energetic 
parameter W = dE/dAo (here M is the oxygen mass) 
and on the characteristic frequency wq of of the 0(2)- 
stretching vibration mode. From the ab-initio-pvo&les 
shown in Fig. MJ^) and Fig. fTUfa). one can derive the 
value of VF = 2.55 eVA-^ which appears to be very close 
to the value 2.4 eVA~^ used in the phenomenological cal- 
culations in Ref. [l!]. Using the value wq ~ 850 cm"^ 
which is in the range of the characteristic frequencies 
of the HSe04 stretching vibrations of Rb3H(Se04)2 re- 
ported in Ref. 37, we obtain the estimate for the pro- 
tonic polaron binding energy: Eq = 0.2 eV. With this 
value, using a high-temperature estimate Ea ~ 5Eq/3 ob- 
tained in the weak proton interaction limit from Refl2TI. 
we find the activation energy for the proton transport 
Ea = 0.33 eV. Table |IT] compares the activation ener- 
gies for the proton conductivity derived on the basis of 
the present DFT-calculations and reported in different 
experimental measurements. The obtained here value 
0.33 eV is very close to the experiment which confirms a 
central role of the 0(2)-distortions for the proton migra- 
tion mechanism. 

It is also worth noting that the snapshot in Fig. H^d) 
corresponds in fact to the two protons located near the 
same ionic group (Se04)^, a configuration similar to the 
high-energy Bjerrum (rotational) defeclj^. The energy 
barrier for the creation of such type of defects is calcu- 
lated from the obtained energy profile and is about 1.5eV 
which is twice as high as the estimated translation barrier 
and which supports negligibly small probability for the 
creation of the Bjerrum defects in this family of proton 
conductors. 



IV. CONCLUSIONS 

We have performed the density functional studies of 
the electronic properties and the mechanisms of hydrogen 
migrations in the system Rb3H(Se04)2 which belongs to 
the technologically promising class of proton conductors 
of M3H(X04)2 crystal family (M=Rb,Cs, NH4; X=S,Se). 
The results of the electronic structure calculations show 
a central role of the lattice dynamics in the process of the 
proton migration. Our findings verify several experimen- 
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tal conclusions obtained in the previous x-ray and NMR 
experiments with Rb3H(Se04)2, which allow to clalify 
the positioning of the vertex oxygens 0(2) during the 
proton migration and to examine the proton positions 
on the hydrogen bond at the intrabond proton transfer. 
The principal conclusion obtained in this work is that 
the free rotational motion of protons without the cor- 
responding complex distortions of Se04 groups cannot 
be considered as an independent part of the Grotthuss 
transport mechanism. In contrast to the almost free pro- 
ton motion implemented in a wide range of the models 
of proton transport, in our work we have considered the 
migration of hydrogen as a cooperative process which in- 
volves a relaxation of the atomic positions and of the 
distances between the different ionic groups XO4. We 
have shown that the migration of the hydrogens is con- 
nected with the significant displacement of the vertex 
oxygens 0(2). These displacements play a decisive role 
in the formation of the protonic polaron and in the con- 
sequent decrease of the activation energy for the proton 
conductivity, a concept introduced in the previous phe- 
nomenological description of the proton transport. The 
present DFT results allow to obtain the first-principle es- 
timates of the hydrogen-bond lengths and energetic pa- 
rameters involved in the existing phenomenological mod- 
eling of the proton transport and to prove the validity of 
the protonic polaron model. Due to their key role in the 
formation of new hydrogen bonds and in the establishing 
the sucessfuU intrabond proton transfer, the oxygen dis- 
placements should be consideted as essential feature of 
the obtained revised mechanism of the hydrogen trans- 
port. The developed approach can be also implemented 
for the modelling of a wide range of the proton conduc- 
tors which in our view substantially increases the degree 
of scientific and possible technological applications of the 
obtained results. 
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